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A conserving T-matrix approximation is applied to the two-dimensional attractive Hubbard model 
in the low-density regime. A set of self-consistent equations is solved in the real- frequency domain 
to avoid the analytic continuation procedure. By tuning the chemical potential the particle density 
was varied in the limits 0.01 < n < 0.18. For the value of the attractive potential U = 8t the 
binding energy of pairs monotonically decreases with increasing n, from its zero-density limit 2.3 1 
and vanishes at a critical density n cr — 0.19. A pairing-induced pseudogap in the single-particle 
density of states is found at low densities and temperatures. 

PACS numbers: 71.10.Fd, 74.20.Mn 



The pseudogap state of underdoped cuprates is cur- 
rently one of the central topics in the studies of high- 
temperature superconductivity (HTSC). The pseudogap 
(PG) is believed to be responsible for marked devia- 
tions from Fermi-liquid behavior observed in a variety of 
probes. Indications from NMR specific heat mea- 

surements H, transport [||,||], optics ||f?]] were followed 
by more recent angle-resolved photoemission Q , Raman 
H and core-electron photoemission pi studies. A con- 
sensus seems to be emerging that the PG is of crucial 
importance for understanding HTSC because it is inti- 
mately related to superconductivity. The PG has the 
same symmetry as the superconducting gap |||)],[ll]], ^ 
decreases with hole doping P,p|Jl2[ and has a universal 
doping dependence when measured in units of the maxi- 
mum critical temperature |nj . 

The PG state gets a simple explanation in phenomeno- 
logical models which assume an effective short-range at- 
tractive interaction between holes doped into a Mott- 
Hubbard insulator Jl^,|l4|]. PG is then associated with 
the binding energy A of hole pairs. Superconductivity 
sets in at lower temperatures T c ~ n pa i rs /m pa i rs [^5| . 
Thus, at small hole densities A and T c are well separated. 
Microscopic origin of such an attraction may be different 
varying from pure electronic [fLsf to pure phononic fl5| l 
mechanisms, see also |L7| . 

In this Letter we address the doping dependence of 
the PG. It has already been discussed by Alexandrov, 
Kabanov, and Mott in |l8[ . They found that the binding 
energy of pairs decreases with doping as a result of the 
screening of the long-range electron-phonon interaction 
and renormalization of the inter-hole coupling. Here we 
will show that PG decreases with doping even in the case 
of short-range inter-hole interaction and without cou- 
pling renormalization. This is a pure many-body effect 
which is expected on the following grounds. If the cou- 
pling exceeds the binding threshold then, at low density, 
the carriers form non-overlapping pairs. As the density 



increases pairs begin to overlap and fcrmions, which con- 
stitute pairs, begin to form a Fermi sea. This process in- 
creases the energy of the system per particle hence raising 
the energy of the two-particle level. Under certain condi- 
tions this level could be pushed up to the single-particle 
continuum. Near the critical density the binding energy 
and the life-time of pairs become small and pairs break 
up. 

In order to verify this idea we study the two- 
dimensional attractive Hubbard model which is the sim- 
plest fermion lattice model with attraction. The model 
is defined by Hamiltonian 



\U\ t t m 

kpq 



written in the standard notation, = —2t(cosk x + 
cosfcj,) is the single-particle bare spectrum, \U\ is the 
coupling strength, TV is the total number of sites in the 
system and chemical potential /i determines the average 
particle density. We will be interested in the low-density 
regime of the model. In this regime one can make use of 
the small gas parameter and select only ladder diagrams 
in the diagrammatic representation of the Bethe-Salpeter 
equation fl9| , |20|| . This procedure leads to the T-matrix 
approximation for the Hubbard model [^lf which we for- 
mulate here in the real-frequency domain. The single- 
particle spectral function A(\h,lj) is determined by the 
real and imaginary parts of the self-energy £' and S" 



A(k,u) = 



-2S"(k,w) 



r k -/i)-£'(k, w )] 2 + [£"(k, w )] 2 ' 



(2) 



The self-energy is expressed, in turn, via the real and 
imaginary parts T' and T" of the T-matrix 
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x T"(q, w + wi) [/f(wi) + / fl (w + wi)] 



(4) 



where /f,b(w) = [exp(/3o;) ± l]" 1 are Fermi- and Bose- 
f unctions respectively, and /3 = 1/T is the inverse abso- 
lute temperature. Finally, the T-matrix is expressed via 
the spectral function as follows 



T(q,w) = 



-1^1 



dwi B(q,wi) .|J7| 

+i — B(q,w) 



(5) 



B(q,w) 
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The integrals with singular kernels in Eqs.([|) and (j5|) 
are understood in the principal value sense. The set of 
equations (^|)-(||) is to be solved self-consistently for given 
|J7|, ji and temperature T, the particle density given by: 



n 
2 



N 4-^ 



du> 
2^ 



(J) 



Single-particle energies are associated with sharp peaks 
in the spectral function A(k, us). Analogously, the en- 
ergies of two-particle bound states can be inferred from 
peaks of the imaginary part of the T-matrix T"'(q, ui) 
which is proportional to the spectral function of the 
two-particle propagator with total quasi-momentum q, 
A 2 {q,u) = -T"(q,Lu)/\U\. The difference between the 
pair energy and the doubled single-particle energies, par- 
ticularly the dependence of this difference on the density 
of particles, is the main object of this work. 

The method of numerical solution of the set of equa- 
tions is well-known: one has to iterate Eqs. (H)-([|) start- 
ing with some initial guess for the spectral function and 
to use the Fast Fourier Transform (FFT) algorithm to 
speed up the calculation of the integrals. We note that 
in previous analogous studies J22|-^| the equations were 
formulated in the imaginary-time domain which required 
an analytic continuation of the iteration results to real 
frequencies. This ill-posed problem of numerical analytic 
continuation leads to ambiguities in the details of recon- 
structed functions. Therefore we prefer to deal directly 
with the real- frequency equations. We have observed 
that the singular structure of the kernels in Eqs. (g)- 
(0) is not an obstacle for the FFT-integration and does 
not slow down the convergence of the iterative process. 
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FIG. 1. Spectral function A(k, w)/(2tt) (solid lines) and 
A2(k,uj) — — T"(k, lj)/\U\ (dashed lines) for three values of 
lattice momentum, k = (0,0), (n/2,ir/2) and (tt,tt) (from 
left to right). Energy is measured with respect to the on-site 
single-particle level i.e. A(k,u>) are shifted by /i and ^(k, ui) 
by 2 fj, from results of self-consistent solution. The position 
of the chemical potential /i = — 6.4 i is shown by the vertical 
long-dashed line. \U\ — 8t, temperature T = 1.0 1 and the 
corresponding particle density is n — 0.056. 



Our calculations were done on a 64 x 64-site lattice and 
on a uniform mesh of 512 points in the frequency interval 
-20 1 < uj < 30 1. Eqs. (§)-(§) were iterated until the 
relative difference between two successive iterations fell 
below 0.01 (usually 20-25 iterations). Note, that no ar- 
tificial damping factor was used in the iteration process. 
The resulting spectral functions were checked against the 
normalization property J j^A(k, uj) = 1 (for other sum 
rules for the Hubbard model see e.g. 



|23J]) which were 
found to be satisfied within 5% accuracy for each k (nor- 
malization was accurate within 1%). The temperature 
was taken T = 1.0 1. 

Figs. |l| and || show results of iterations for two different 
values of the chemical potential \x = — 6 At and /i = 
— 5.8 i corresponding to particle densities n — 0.056 and 
n = 0.111 respectively. To obtain the binding energy 
of a bound state with total quasi-momentum q one has 
to compare the position of the corresponding peak in 
A 2 (q_, uj) with the doubled position of a peak in spectral 
function A(k, u>) taken at half momentum (q/2). Thus, 
in Fig. |] the binding energy of the pair with q = (0, 0) is 
2J5i(0,0) - E 2 (0,0) w 2.3 1, the one of the pair with q = 
(7r,7r) is 2 £1(71-/2,71-/2) -£ 2 (7r,7r) «8.0t. An analysis of 
Figs, [l] and || and analogous data for seven other densities 
in the region 0.01 < n < 0.18 (not shown) reveals the 
following common trends, i) The single-particle peak as 
a function of momentum k spans an energy interval of 
width k 8t for all the densities studied. This implies 
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FIG. 2. The same as Fig. 1 but for chemical potential 
fi = — 5.8t corresponding to a larger density n = 0.111. 

small effect of the interaction on the shape and width of 
the single-particle energy spectrum. The same conclusion 
has been reached in Ref. |^3| . ii) The lower boundary of 
the single-particle band £1(0,0) shifts towards negative 
values as density increases. However, the shift observed is 
smaller than the mean-field value —\U\n/2. Apparently, 
this is due to the correlation effects. At \U\ = 8t the 
effective radius of a pair is « 1 — 2 lattice spacings, see 
e.g. pfjj . Therefore the interaction of a particle with 
spin a excited to the continuum spectrum, with a bound 
particle of spin — a is effectively screened by the latter's 
spin a partner, iii) The energy of the two-particle state 
increases with the density, as indeed was expected on the 
grounds of the qualitative considerations, iv) The single 
particle and the pair spectral functions become narrower 
at large lattice momenta. The life-time of the bound 
states with large momenta increases which results from 
reduced phase space of the scattering process into two 
single-particle states, v) The binding energy of the pair 
with the highest momentum q = (tt, n) is always close 
to \U\ (which is equal to 8 1 in our case). This result of 
a fully self-consistent treatment is in perfect agreement 
with the exact solution of the two-body problem (the 
zero-density limit) |^7j and the general theory of the 77- 
resonance |2f| according to which a collective mode at 
q = (tt, tt) with binding energy \U\ always exists in the 
attractive Hubbard model. 

At low density only states with small momenta are 
occupied. Therefore we now concentrate on the prop- 
erties of pairs with q = (0,0). We have analysed the 
doping-dependence of the binding energy using solutions 
of the self-consistent equations for nine different densi- 
ties in the interval 0.01 < n < 0.18. For the smallest 
density studied, n = 0.01, the binding energy was found 
to be w 2.3 i. For the two-particle problem the binding 
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FIG. 3. Pair binding energy vs. particle density. The 
points are determined from the position of the maxima of 
peaks of spectral functions. Bars indicate the width of peaks. 

energy in the ground state is A(n = 0) = 2.1 1 ^?J. The 
proximity of these two values confirms the correctness of 
our numerical treatment of the self-consistent equations. 
The small difference is apparently due to the finite-size 
effects. The density dependence of the binding energy is 
shown in Fig.||. It is a monotonically decreasing func- 
tion which vanishes at n cr w 0.19. The points have been 
calculated from the positions of the maxima of one- and 
two-particle spectral functions. The peaks' widths are 
associated with finite life-time of the excitations. The 
bars in Fig.|| indicate the width (the inverse life-time) of 
the two-particle spectral weight. Obviously, if the energy 
gap between the bound states and the single-particle con- 
tinuum diminishes the scattering process of the pair into 
two single-particle states intensifies. This should lead to 
decrease of the life-time and widening of the peaks. This 
is exactly what is seen in the self-consistent results. We 
have checked that the decrease of the binding energy is 
to be attributed mainly to the increase of the energy of 
two-particle states. There is also a contribution from the 
negative shift of single-particle continuum, but we have 
found this process to account for « 10% of the whole 
effect. It is also clear that the region of stability of the 
bound states depends on the coupling constant: it ex- 
pands for larger |U|'s and shrinks for smaller ones. We 
have repeated all the calculations for \U\ = 4.0 1 and 
T = 0.2 1 and found that for density n — 0.10 the pair 
with q = (0, 0) is not a well-defined excitation, its bind- 
ing energy and life-time being too small. The pair ex- 
citation separates from the single-particle continuum at 
q w (7r/2,7r/2). At the corner of the Brillouin zone the 
pair's binding energy is ss \U\ = 4.0 1. Therefore, the cou- 
pling is to be increased if one requires well-defined pair 
excitations at q = (0,0). At large couplings the conver- 
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FIG. 4. The low-energy part of the density of states 
for \U\ — 8t. Solid lines are for T = 0.5 t and densities 
n = 1.0,1.7,3.3 and 6.3% (from the lowest curve upwards). 
Dashed lines are n — 1.7% and T — 0.4 i (upper curve) and 
T — 0.3 t. Vertical lines mark the energies coo such that 
X^l dup{u) = n. 

gence of the iterative process worsens but restores with 
increasing temperature. Thus, we were able to trace the 
destruction of pairs in full at \U\ — 8.0 1 and T = 1.0 1. 

The density dependence of the pair binding energy 
found in our calculations is similar to that of the pseudo- 
gap observed in the underdoped regime of high-T c super- 
conductors. This supports the real-space pairing as the 
cause for the pseudogap phenomenon. In order to test 
and clarify this relation further we calculated the single- 
particle density of states p(ui) from the solution of the T- 
matrix equations. The low-energy part of p(ui) is shown 
m Fig. |. At small n we observe an extra peak in p{u>) 
which is apparently due to the presence of bound states, 
since its separation from the main part of p(u>) is approx- 
imately the binding energy. (Compare with the atomic 
limit of the Hubbard model, see e.g. ^3|.) At larger n 
the feature shifts to higher energies and finally disappears 
thereby qualitatively repeating the n-dependence of the 
binding energy discussed above. Two comments are to 
be added. First, the density of states is an integral of the 
single-particle spectral function over momenta k. Since 
at different k the additional maxima appear at different 
energies the feature is less pronounced in p{us) than in 
individual A(k, ui). Second, upon lowering the tempera- 
ture the local minimum in p(u>) deepens, see the dashed 
curves in Fig. |j. This implies the opening of a pseudogap 
at low temperatures, as experimentally observed in the 
underdoped cuprates. Note that the position of the peak 
is not affected by temperature which confirms that it is 
the pair binding energy. 

In conclusion, we have studied the low-density regime 



of the two-dimensional attractive Hubbard model by 
means of the self-consistent T-matrix approximation. 
The set of self-consistent equations was solved in the real- 
frequency domain to avoid uncertainties related to the 
analytic continuation procedure. The binding energy of 
pairs was found to be a monotonically decreasing func- 
tion of the total density. We interprete this as a pure 
many-body effect resulting from the rise of two-particle 
levels to the single-particle continuum due to the pack- 
ing effect. We believe this effect is generic to low-density 
regimes of fermionic models with attraction and it does 
not depend on the dimensionality or particular details of 
the interaction. In particular it should remain in more 
realistic models of low-density hole systems such as the 
nearest-neighbors attraction model with ci-symmetrical 
bound states. We observe also a pairing-induced pseudo- 
gap in the density of states with the characteristic size 
equal to the binding energy of pairs. The pseudogap de- 
velops at low densities upon lowering temperature. 
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